In [1]:
%matplotlib inline
%config InlineBackend.figure_format = 'retina'

import numpy as np
import pandas as pd
import seaborn as sns
import scanpy as sc
from glob import iglob
import scrublet as scr
import anndata
import scanpy.external as sce
import scvelo as scv
#import scrublet as scr # requires 'pip install scrublet'
import os
import sklearn
import bbknn
from sklearn.linear_model import LogisticRegression
import matplotlib as mpl
import scipy
import matplotlib.pyplot as plt
import pickle

np.random.seed(0)
sc.settings.verbosity = 3  # verbosity: errors (0), warnings (1), info (2), hints (3)
sc.logging.print_versions()
sc.settings.set_figure_params(dpi=100,dpi_save=800)  # low dpi (dots per inch) yields small inline figures
-----
anndata     0.8.0
scanpy      1.9.1
-----
OpenSSL                                     21.0.0
PIL                                         9.0.1
aa8f2297d25b4dc6fd3d98411eb3ba53823c4f42    NA
absl                                        NA
annoy                                       NA
asciitree                                   NA
asttokens                                   NA
astunparse                                  1.6.3
backcall                                    0.2.0
bbknn                                       NA
bcrypt                                      3.2.0
beta_ufunc                                  NA
binom_ufunc                                 NA
boto3                                       1.21.32
botocore                                    1.24.32
bottleneck                                  1.3.4
brotli                                      NA
cairo                                       1.21.0
certifi                                     2021.10.08
cffi                                        1.15.0
chardet                                     4.0.0
charset_normalizer                          2.0.4
cloudpickle                                 2.0.0
colorama                                    0.4.4
cryptography                                3.4.8
cycler                                      0.10.0
cython_runtime                              NA
cytoolz                                     0.11.0
dask                                        2022.02.1
dateutil                                    2.8.2
debugpy                                     1.5.1
decorator                                   5.1.1
defusedxml                                  0.7.1
dot_parser                                  NA
entrypoints                                 0.4
executing                                   0.8.3
fasteners                                   0.17.3
flatbuffers                                 NA
fsspec                                      2022.02.0
gast                                        NA
google                                      NA
h5py                                        3.6.0
idna                                        3.3
igraph                                      0.9.11
ipykernel                                   6.9.1
ipython_genutils                            0.2.0
ipywidgets                                  7.6.5
jedi                                        0.17.2
jinja2                                      3.1.2
jmespath                                    0.10.0
joblib                                      1.1.0
jupyter_server                              1.13.5
keras                                       2.9.0
kiwisolver                                  1.3.2
leidenalg                                   0.8.10
llvmlite                                    0.38.0
louvain                                     0.7.1
markupsafe                                  2.0.1
matplotlib                                  3.5.1
matplotlib_inline                           NA
mkl                                         2.4.0
mpl_toolkits                                NA
msgpack                                     1.0.2
natsort                                     8.1.0
nbinom_ufunc                                NA
nt                                          NA
ntsecuritycon                               NA
numba                                       0.55.1
numcodecs                                   0.10.0
numexpr                                     2.8.1
numpy                                       1.21.5
opt_einsum                                  v3.3.0
packaging                                   21.3
pandas                                      1.4.2
parso                                       0.7.1
pickleshare                                 0.7.5
pkg_resources                               NA
prompt_toolkit                              3.0.20
psutil                                      5.8.0
pure_eval                                   0.2.2
pycparser                                   2.21
pydev_ipython                               NA
pydevconsole                                NA
pydevd                                      2.6.0
pydevd_concurrency_analyser                 NA
pydevd_file_utils                           NA
pydevd_plugins                              NA
pydevd_tracing                              NA
pydot                                       1.4.2
pygments                                    2.11.2
pynndescent                                 0.5.7
pyparsing                                   3.0.4
pythoncom                                   NA
pytz                                        2021.3
pywintypes                                  NA
requests                                    2.27.1
scipy                                       1.7.3
scrublet                                    NA
scvelo                                      0.2.4
seaborn                                     0.11.2
session_info                                1.0.0
setuptools                                  61.2.0
setuptools_scm                              NA
six                                         1.16.0
sklearn                                     1.0.2
snappy                                      NA
socks                                       1.7.1
sphinxcontrib                               NA
stack_data                                  0.2.0
statsmodels                                 0.13.2
tblib                                       1.7.0
tensorboard                                 2.9.1
tensorflow                                  2.9.1
termcolor                                   1.1.0
texttable                                   1.6.4
threadpoolctl                               2.2.0
tlz                                         0.11.0
toolz                                       0.11.2
tornado                                     6.1
tqdm                                        4.64.0
traitlets                                   5.1.1
typing_extensions                           NA
umap                                        0.5.3
urllib3                                     1.26.9
wcwidth                                     0.2.5
win32api                                    NA
win32com                                    NA
win32security                               NA
wrapt                                       1.12.1
yaml                                        6.0
zarr                                        2.12.0
zipp                                        NA
zmq                                         22.3.0
zope                                        NA
-----
IPython             8.2.0
jupyter_client      6.1.12
jupyter_core        4.9.2
jupyterlab          3.3.2
notebook            6.4.8
-----
Python 3.9.12 (main, Apr  4 2022, 05:22:27) [MSC v.1916 64 bit (AMD64)]
Windows-10-10.0.19044-SP0
-----
Session information updated at 2022-09-07 09:01
In [2]:
adata = sc.read_h5ad('./01.single data/scRNA_qc.h5ad')
In [3]:
sc.pp.normalize_total(adata, target_sum=1e4)
sc.pp.log1p(adata)
normalizing counts per cell
    finished (0:00:00)
In [4]:
sc.pp.highly_variable_genes(adata, min_mean=0.0125, max_mean=3, min_disp=0.5)
print("Highly variable genes: %d"%sum(adata.var.highly_variable))
sc.pl.highly_variable_genes(adata)
extracting highly variable genes
    finished (0:00:01)
--> added
    'highly_variable', boolean vector (adata.var)
    'means', float vector (adata.var)
    'dispersions', float vector (adata.var)
    'dispersions_norm', float vector (adata.var)
Highly variable genes: 2987
In [5]:
sc.pp.regress_out(adata, ['n_counts', 'pct_counts_mt'])
sc.pp.scale(adata, max_value=10)
regressing out ['n_counts', 'pct_counts_mt']
    sparse input is densified and may lead to high memory use
    finished (0:22:29)
In [6]:
sc.tl.pca(adata, svd_solver='arpack')
sc.pl.pca_variance_ratio(adata, log=True,n_pcs = 50)
sns.lineplot(data = adata.uns['pca']['variance_ratio'],palette= 'tab10',linewidth = 3)
adata.uns['pca']['variance']
sns.lineplot(data =adata.uns['pca']['variance'],palette='tab10',linewidth = 3)
computing PCA
    on highly variable genes
    with n_comps=50
    finished (0:00:11)
Out[6]:
<AxesSubplot:>
In [7]:
%%time
# trying out harmony here
# Check PCA again - needed for harmony - and extract PCA matrix and batch array
n_pcs = 30
pca = adata.obsm['X_pca'][:, 0:n_pcs]
batch = adata.obs['ID']
CPU times: total: 0 ns
Wall time: 0 ns
In [8]:
# Batch-correct the PCA using HARMONY method
%load_ext rpy2.ipython
In [9]:
%%time
%%R -i pca -i batch -o hem

library(harmony)
library(magrittr)

hem <- HarmonyMatrix(pca, batch, theta=1, do_pca=FALSE)
hem = data.frame(hem)
Exception ignored from cffi callback <function _consolewrite_ex at 0x000001FB01BE33A0>:
Traceback (most recent call last):
  File "D:\conda\lib\site-packages\rpy2\rinterface_lib\callbacks.py", line 133, in _consolewrite_ex
    s = conversion._cchar_to_str_with_maxlen(buf, n, _CCHAR_ENCODING)
  File "D:\conda\lib\site-packages\rpy2\rinterface_lib\conversion.py", line 138, in _cchar_to_str_with_maxlen
    s = ffi.string(c, maxlen).decode(encoding)
UnicodeDecodeError: 'utf-8' codec can't decode byte 0xd4 in position 0: invalid continuation byte
R[write to console]: Harmony 1/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony 2/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony 3/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony 4/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony 5/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony 6/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony 7/10

R[write to console]: 0%   10   20   30   40   50   60   70   80   90   100%

R[write to console]: [----|----|----|----|----|----|----|----|----|----|

R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: *
R[write to console]: |

R[write to console]: Harmony converged after 7 iterations

CPU times: total: 4min 51s
Wall time: 4min 51s
In [10]:
# Add harmony values to the anndata object
adata.obsm['X_pca'] = hem.values
In [16]:
for i in range(10,60,5):
    sc.pp.neighbors(adata,n_pcs=30, n_neighbors= i)
    sc.tl.umap(adata)
    print(f'neighbors的数值为:{i}')
    sc.pl.umap(adata)
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:02)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:00:43)
neighbors的数值为:10
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:03)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:00:43)
neighbors的数值为:15
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:03)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:00:47)
neighbors的数值为:20
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:04)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:00:51)
neighbors的数值为:25
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:05)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:00:55)
neighbors的数值为:30
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:06)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:00:57)
neighbors的数值为:35
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:07)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:01:00)
neighbors的数值为:40
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:09)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:01:02)
neighbors的数值为:45
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:11)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:01:04)
neighbors的数值为:50
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:14)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:01:07)
neighbors的数值为:55
In [15]:
sc.pp.neighbors(adata,n_pcs=30, n_neighbors= 55)
sc.tl.umap(adata)
sc.pl.umap(adata)
computing neighbors
    using 'X_pca' with n_pcs = 30
    finished: added to `.uns['neighbors']`
    `.obsp['distances']`, distances for each pair of neighbors
    `.obsp['connectivities']`, weighted adjacency matrix (0:00:14)
computing UMAP
    finished: added
    'X_umap', UMAP coordinates (adata.obsm) (0:01:06)
In [17]:
adata1 = sc.read_h5ad('./01.single data/scRNA_celltype.h5ad')
In [23]:
umap = adata1.obsm['X_umap']
In [24]:
adata.obsm['X_umap'] = umap
In [25]:
sc.pl.umap(adata)
In [26]:
for i in range(1,16):
    sc.tl.leiden(adata,resolution=round(i*0.1,1),key_added=f'cluster_{round(i*0.1,1)}')
running Leiden clustering
    finished: found 13 clusters and added
    'cluster_0.1', the cluster labels (adata.obs, categorical) (0:02:05)
running Leiden clustering
    finished: found 17 clusters and added
    'cluster_0.2', the cluster labels (adata.obs, categorical) (0:01:45)
running Leiden clustering
    finished: found 21 clusters and added
    'cluster_0.3', the cluster labels (adata.obs, categorical) (0:04:08)
running Leiden clustering
    finished: found 22 clusters and added
    'cluster_0.4', the cluster labels (adata.obs, categorical) (0:03:00)
running Leiden clustering
    finished: found 22 clusters and added
    'cluster_0.5', the cluster labels (adata.obs, categorical) (0:03:35)
running Leiden clustering
    finished: found 24 clusters and added
    'cluster_0.6', the cluster labels (adata.obs, categorical) (0:05:13)
running Leiden clustering
    finished: found 25 clusters and added
    'cluster_0.7', the cluster labels (adata.obs, categorical) (0:03:48)
running Leiden clustering
    finished: found 25 clusters and added
    'cluster_0.8', the cluster labels (adata.obs, categorical) (0:04:42)
running Leiden clustering
    finished: found 26 clusters and added
    'cluster_0.9', the cluster labels (adata.obs, categorical) (0:03:26)
running Leiden clustering
    finished: found 28 clusters and added
    'cluster_1.0', the cluster labels (adata.obs, categorical) (0:02:34)
running Leiden clustering
    finished: found 30 clusters and added
    'cluster_1.1', the cluster labels (adata.obs, categorical) (0:03:19)
running Leiden clustering
    finished: found 31 clusters and added
    'cluster_1.2', the cluster labels (adata.obs, categorical) (0:04:36)
running Leiden clustering
    finished: found 33 clusters and added
    'cluster_1.3', the cluster labels (adata.obs, categorical) (0:05:04)
running Leiden clustering
    finished: found 33 clusters and added
    'cluster_1.4', the cluster labels (adata.obs, categorical) (0:08:35)
running Leiden clustering
    finished: found 36 clusters and added
    'cluster_1.5', the cluster labels (adata.obs, categorical) (0:04:41)
In [27]:
sc.pl.umap(adata,color = [f'cluster_{round(i*0.1,1)}' for i in range(1,16)],ncols = 3,wspace=0.5)
In [28]:
sc.pl.umap(adata, color=['Type','cluster_0.6'],add_outline=True, legend_loc='right margin',size=15,
               legend_fontsize=15, legend_fontoutline= 1,frameon= True,wspace=0.3)
In [29]:
adata.obs['Type'].value_counts()
Out[29]:
Control    45618
Vehicle    16099
Drug       13528
Taxol      12831
Name: Type, dtype: int64
In [21]:
adata.obs['Type'].value_counts()
Out[21]:
Control    45618
Vehicle    16099
Drug       13528
Taxol      12831
Name: Type, dtype: int64
In [30]:
sc.tl.rank_genes_groups(adata,'cluster_0.6',method='wilcoxon')
sc.pl.rank_genes_groups(adata, n_genes=25, sharey=False)
ranking genes
    finished: added to `.uns['rank_genes_groups']`
    'names', sorted np.recarray to be indexed by group ids
    'scores', sorted np.recarray to be indexed by group ids
    'logfoldchanges', sorted np.recarray to be indexed by group ids
    'pvals', sorted np.recarray to be indexed by group ids
    'pvals_adj', sorted np.recarray to be indexed by group ids (0:12:11)
In [31]:
pd.DataFrame(adata.uns['rank_genes_groups']['names']).head(35)
Out[31]:
0 1 2 3 4 5 6 7 8 9 ... 14 15 16 17 18 19 20 21 22 23
0 Apoe Cst3 Gramd3 Plp1 Cd74 S100a9 Ctla4 Srgn Dusp9 Plcb1 ... Krt71 Emb Col1a2 Plp1 Syne1 Lsamp Siglech Hbb-bs Sptbn1 Itgb8
1 Ctsd Ctss Ptpn22 Mbp H2-Aa S100a8 Icos Clec4d Gm37035 S100a4 ... Gm47501 Rora Sparc Apod Fas Opcml Bst2 Hba-a1 Flt1 Timp3
2 C1qa Basp1 Itk Pcdh9 H2-Ab1 Lcn2 Rora A630073D07Rik Gm44780 Gda ... Gm10610 Ikzf3 Col1a1 Mbp Dbi Dscam St8sia4 Hbb-bt Pecam1 Utrn
3 C1qb C1qb Cblb Nkain2 H2-Eb1 Ngp Ets1 Vmn2r117 Gm13855 Ifitm3 ... Gm48923 Icos Dcn Gm42418 Klf5 Lhfpl3 Tcf4 Hba-a2 Epas1 Stard13
4 Trem2 Zeb2 Vps37b Mal Lyz2 Wfdc21 Tnfrsf4 Gm2244 Gm15792 F13a1 ... Gm4858 Zbtb16 Gsn Pcdh9 Cd24a Nrxn1 Aff3 Gypa Ptprb Dmd
5 C1qc Fnip2 Ets1 Aplp1 Lgals3 Pglyrp1 Tnfrsf9 Nat1 Vmn2r4 Plac8 ... 4930527G23Rik Camk4 Mgp Nkain2 Gm47708 Vcan Ighm Ermap Adgrl4 Nkain2
6 Ctss Rbm47 Grap2 Trf Tgfbi Ifitm6 Odc1 4930533N22Rik Gm42697 Emilin2 ... Gm40638 Il7r Bgn Ptgds Gm15082 Tnr Tspan13 Hmbs Igfbp7 Sparc
7 Hexb Atf3 Mbnl1 Cnp Cybb Syne1 Kdm2b Ighv3-6 Olfr1152 Lyz2 ... Cpn2 Dgat1 Cfh Camk1d Klhl1 Nlgn1 Plac8 Mki67 Nfib Mpz
8 Cst3 Csf1r Trbc2 Scd2 Psap Hdc Il2rb Gm29674 9230009I02Rik Ms4a6c ... Ang2 St6galnac3 Cped1 Cldn11 Gm48821 Lrrc4c Mef2c Ank1 Abcg2 Lmna
9 Cd14 Lrmda Camk4 Cldn11 Sdc4 Camp Camk4 5033428I22Rik Dkk1 Thbs1 ... Gm26654 Cxcr6 Serpinh1 Car2 Vmn1r38 Cadm2 Rell1 Slc4a1 Mecom Dag1
10 Ctsz Lgmn Zc3hav1 Magi2 Mafb Anxa1 Itk AL845275.1 Vmn2r83 Tgfbi ... Igkv3-5 Zeb1 Igfbp7 Scd2 Gm36525 Pcdh9 Lrp8 Car2 Ptprg Hspa12a
11 Ctsb C1qc Skap1 Car2 Cstb G0s2 Hif1a Gm11636 Gm11924 Alox5ap ... Gm14635 Tex2 Col3a1 Mobp Gm48155 Pcdh7 Ccr9 Tfrc Utrn Sorcs1
12 Timp2 Epb41l2 Prkca Slc24a2 Ctsc Sgms2 Crem Gm35281 Lypd8 Vim ... Gm20268 Crem Id3 Mal Gm34350 Adgrl3 Eepd1 Prdx2 Bsg Plekha4
13 Spp1 C1qa Gimap6 Mag Pid1 Retnlg Inpp4b Gm47486 4930470H14Rik Ly6c2 ... Ddi1 Igf1r Serpinf1 Aplp1 Wfdc15a Sox6 Mctp2 Kel Cdh5 Plp1
14 Cd63 Hexb Gm2682 Apod Ms4a6c Cd177 Tnfrsf18 Gm15856 Slc51b Fn1 ... Igkv14-111 D16Ertd472e Nedd4 Trf Gm44204 Cntn1 Ly6d Slc25a21 Egfl7 Ahnak
15 Tyrobp Ctsb Pde7a Qk Cxcl16 Msra Akap13 Slc12a1 Gm31138 Lgals3 ... Ighv8-10 Itk Fstl1 Tmeff2 2610318M16Rik Ptprz1 Gpr171 Sox6 Eng Qk
16 Abca1 Ctsd Ifngr1 Mobp Fcgr2b S100a11 Trac Trim80 4930402F11Rik Vav3 ... Npbwr1 Nr3c1 Itm2a Magi2 Trav13d-1 Grid2 Lgals1 Cd24a Tcf4 Afap1l2
17 Cd68 Tanc2 Prkch Edil3 Pla2g7 Ltf S100a10 Gm49497 Gm31786 Ccl9 ... Fam162b Vps37b Zbtb20 S100a9 Gm50114 Adgrb3 Tex2 Hemgn Cldn5 Sorbs1
18 Ctsl Tns3 Kdm2b Cryab Anxa5 Rasa2 Cd2 Eras Gm13206 F10 ... Gm16144 Il18r1 Cald1 Slc24a2 4930556G01Rik Plcb1 Sell Slc25a37 Erg Megf9
19 Fcrls Mafb Gimap4 Sept4 Mmp14 Hp Ndfip1 Fam43b Gm29825 Arhgap26 ... 1700001G01Rik Il2rb Gpc6 Mag 4930458B22Rik Erbb4 Blnk Fech Itm2a Hspg2
20 Syngr1 Pmp22 Zeb1 Dbndd2 Ftl1 Plaur Rabgap1l Gm40038 Lipo4 S100a6 ... Gm26650 Nrip1 Col5a2 Cnp Adam28 Pcdh15 Cd7 Ctse Ly6c1 Lama2
21 Ftl1 C3ar1 Ms4a4b Tmeff2 Npc2 Adpgk Cd3g Gm38661 Gm11642 Ccr2 ... Gm49710 Kdm2b Igfbp6 AY036118 Gm17266 Ppfibp1 Selplg Cpox Ebf1 Tead1
22 Csf1r Apoe Emb Plcl1 Vim Arhgdib Cd3d Gm16143 AB041806 Ifitm2 ... Kcna10 Ramp3 Mmp14 Edil3 A930002I21Rik Mdga2 Irf8 Spta1 Tjp1 Ank3
23 Dock4 Lhfpl2 Ablim1 Ermn Sdcbp Mcemp1 Itpr1 Gm5458 Gm15245 Napsa ... Gm28410 Ets1 Serping1 Apoc1 Gm49124 Nav3 Ly6e Sptb Ppfibp1 Sptbn1
24 Mpeg1 Cadm1 Cd3d Lrp1b Prdx1 Gsr Dgat1 Gm11264 Ighv9-3 Ahnak ... Gm49932 Il23r Bicc1 Ptprd Gm10646 Lrp1b Cytip Rhd Wwtr1 Auts2
25 Cd83 Trem2 Cd28 Ptprd Ninj1 Sp140 Aebp2 Dcaf12l2 Gm14435 Cd44 ... Krtap28-13 Aebp2 Egr1 Qrfprl Gm45238 Npas3 Ncf1 Aqp1 Arhgap29 Art3
26 Fth1 Hpgds Tmsb10 Enpp2 Ctss Itgam Trbc2 Chrdl2 Defb43 Cybb ... Gm35147 Tmem64 Ahnak Gm17206 A330094K24Rik Anks1b Glcci1 Alas2 Tshz2 Dst
27 Grn Ptprm Fryl Kif1b Fth1 Antxr2 H2-K1 Gm14664 Gm15179 S100a10 ... Avpr1b Cd164 Cebpd Gm13209 Gkn3 Il1rapl1 Flt3 Gpx1 Cyyr1 Kcna1
28 Lgmn Tgfbr1 Tut4 Gpm6b Ctsb Dgat1 Tmsb10 Zfp92 Spz1 Esd ... Gm35974 Arap2 Pcolce Gm36107 Traj27 Csmd1 Sub1 Rhag Adgrf5 Tenm3
29 Fcer1g Blnk H2-Q7 Mast4 H2-DMa Igf1r Cd28 Trpv1 Ccnb3 Actg1 ... Gm49156 Fryl Klf9 9530036O11Rik 1700010B13Rik Rtn1 Cox6a2 Abcb10 Itga6 Pmp22
30 Ctsa Rab7b Dock10 Syt11 Cxcl2 Alox5ap Cytip Akr1cl Kbtbd13 Vcan ... Gm36231 Ccr2 Sptbn1 Glra3 Gm15668 Gria4 Snx9 Metap2 Zbtb20 Dlc1
31 Bcl2a1b Ccl4 Gimap3 Pde4b Ahnak Lmnb1 Cd4 Gm33125 Gm15891 Gsr ... Foxe1 Ppp1cc Fgfr1 Uncx Ube2d2b Ncam2 Rabgap1l E2f2 Timp3 Lamb1
32 Lyz2 Cd83 Stat4 St18 Ms4a6d Tbc1d8 Prkca Gm14565 Cyp2c66 Msrb1 ... Gm30409 Ifngr1 Pbx1 Gm11682 Nkx2-4 Brinp3 Ifnar2 Snca Fry Vcl
33 Itgb5 Ly86 Stk17b Tspan2 Cd14 Prdx5 S100a4 4930568D16Rik Gm48306 Trps1 ... Gm39329 Ern1 Col6a1 9530026F06Rik Trav9-1 Nxph1 Sla2 Tspan33 Cd200 Gpm6b
34 Bcl2a1a Snx24 Rundc3b Sec11c Fnip2 Mxd1 D16Ertd472e Gm48703 Gm16180 Clec4d ... Gm28867 Prelid2 Rcn3 Gm6588 Olfr729 Luzp2 Snx5 Nxpe2 Ly6a Serinc5

35 rows × 24 columns

In [32]:
result = adata.uns['rank_genes_groups']
groups = result['names'].dtype.names
res = pd.DataFrame(    {group + '_' + key: result[key][group]    for group in groups for key in ['names', 'pvals','logfoldchanges','pvals_adj','scores']})
In [33]:
res.to_csv("diff_cluster.csv")
In [36]:
sc.pl.umap(adata, color=['cluster_0.6'],add_outline=True, legend_loc='on data',size=18,title='',
               legend_fontsize=10, legend_fontoutline= 1,frameon= True)
In [281]:
markers = ["Apoe","Ctss",#Microglial
           "Cd3d",#Lymphocyte
           "Rora",#Th17
           "Nkg7",#NK
           "Cd79a","Cd79b",#B
           "Lyz2","Cd68",#Macrophage
           "S100a8","S100a9",#Neutrophil
           'Ccr7',"Cd74",#Dendritic cell
           "Ly6c2","Csf1r",#Monocyte
           "Aqp4","Gfap",#Astrocyte
           "Dcn","Lum",#Fibroblast
           "Map2","Dcx",#Neuron
           "Tmem212","Foxj1",#Ependymal cell
           "Plp1","Opalin",#Oligodendrocyte
           "Olig2",#Oligodendrocyte precursor cell少突胶质前体细胞
           'Cldn5',"Pecam1",#endothelial cell
           "Dcn","Lum",#Fibroblast
           "Cdk1","Mki67",#Dividing cell
          ]
sc.pl.dotplot(adata,markers,'cluster_0.6',cmap='OrRd',swap_axes=True)
In [267]:
sc.settings.set_figure_params(dpi=100,dpi_save=600,figsize=(5,5))
sc.pl.umap(adata, color=['cluster_0.6','Type'],add_outline=True, legend_loc='on data',size=18,
               legend_fontsize=10, legend_fontoutline= 1,frameon= True,wspace=0.3)
In [279]:
sc.pl.umap(adata,color=["S100a9","Gfap","C1qc","Cd14"],ncols=4,cmap='OrRd',size=15,legend_fontsize=12)
In [358]:
def annot_prelim(cluster):
    if cluster in ["0","2","8"]:
         return('Microglial')
    elif cluster in ['1','7','14']:
        return('Lymphocyte')
    elif cluster in ['3']:
        return('OPC')
    elif cluster in ['4','9']:
        return('Macrophages')
    elif cluster in ['5','6','15']:
        return('Neutrophil')
    elif cluster =='10':
        return('Dendritic cell')
    elif cluster =='11':
        return('NK cell')
    elif cluster =='12':
        return('B cell')
    elif cluster =='13':
        return('Monocyte')
    elif cluster =='16':
        return('Fibroblast')
    elif cluster =='17':
        return('Astrocyte')
    elif cluster =='18':
        return('ODC')
    elif cluster =='19':
        return('Ependymal')
    elif cluster =='20':
        return('Neuron')
    elif cluster =='21':
        return('Dividing cell')
    elif cluster =='22':
        return('Endothelial')
    else:
        return('Th17')
In [359]:
adata.obs['Celltype'] = adata.obs['cluster_0.6'].map(annot_prelim)
In [360]:
sc.settings.set_figure_params(dpi=100,dpi_save=600,figsize=(6,6))
sc.pl.umap(adata, color='Celltype',add_outline=True, legend_loc='right margin',
               legend_fontsize=15, legend_fontoutline= 1,frameon= True)
In [361]:
adata.obs["Celltype"] = adata.obs.Celltype.cat.set_categories([
       'Microglial','Astrocyte', 'OPC','ODC',
       'Neuron','Ependymal','Fibroblast','Endothelial',
       'Dividing cell','Lymphocyte','NK cell','Th17',
       'B cell','Monocyte','Macrophages','Neutrophil',
       'Dendritic cell'])
In [362]:
adata.obs["Celltype"] = adata.obs["Celltype"].astype('category')
adata.obs["Celltype"].cat.categories
Out[362]:
Index(['Microglial', 'Astrocyte', 'OPC', 'ODC', 'Neuron', 'Ependymal',
       'Fibroblast', 'Endothelial', 'Dividing cell', 'Lymphocyte', 'NK cell',
       'Th17', 'B cell', 'Monocyte', 'Macrophages', 'Neutrophil',
       'Dendritic cell'],
      dtype='object')
In [363]:
new_colors = np.empty(len(adata.obs["Celltype"].cat.categories), dtype=object) 
new_colors[[0]] = '#FF0033' # Microglial
new_colors[[1]] = '#339933' # Astrocyte
new_colors[[2]] = '#FF6600' # 'OPC'
new_colors[[3]] = '#FF9999' # Oligodendrocyte
new_colors[[4]] = '#990066' # Neuron
new_colors[[5]] = '#CC00CC' #Ependymal
new_colors[[6]] = '#33FF33' # Fibroblast
new_colors[[7]] = '#000099' # Endothelial
new_colors[[8]] = '#FFFF00' # Dividing cell
new_colors[[9]] = '#FF6666' # Lymphocyte
new_colors[[10]] ='#FFCC66' # NK cell
new_colors[[11]] ='#990000' # Th17
new_colors[[12]] ='#FF9933' # B cell
new_colors[[13]] ='#0099CC' # Monocyte
new_colors[[14]] ='#0066CC' # Macrophages
new_colors[[15]] ='#9999FF' # Neutrophil
new_colors[[16]] ='#3399FF' # Dendritic cell
adata.uns["Celltype_colors"] = new_colors
new_colors
Out[363]:
array(['#FF0033', '#339933', '#FF6600', '#FF9999', '#990066', '#CC00CC',
       '#33FF33', '#000099', '#FFFF00', '#FF6666', '#FFCC66', '#990000',
       '#FF9933', '#0099CC', '#0066CC', '#9999FF', '#3399FF'],
      dtype=object)
In [364]:
sc.settings.figdir = './result/'
sc.settings.set_figure_params(dpi=100,dpi_save=600,figsize=(8,8))
sc.pl.umap(adata, color='Celltype',legend_loc='on data',size=5,
           legend_fontsize=8, legend_fontoutline= 2,title='Cell Type',add_outline=True,
              frameon= True,save='celltype.pdf')
WARNING: saving figure to file result\umapcelltype.pdf
In [365]:
%load_ext rpy2.ipython
The rpy2.ipython extension is already loaded. To reload it, use:
  %reload_ext rpy2.ipython
In [366]:
%%R 
library(future)
library(RColorBrewer)
library(viridis)
library(Seurat)
library(harmony)
library(tidyverse)
library(patchwork)
library(ggplot2)
library(cowplot)
library(dplyr)
library(data.table)
In [367]:
adata
Out[367]:
AnnData object with n_obs × n_vars = 88076 × 25375
    obs: 'ID', 'Type', 'batch', 'n_genes', 'n_genes_by_counts', 'total_counts', 'total_counts_mt', 'pct_counts_mt', 'total_counts_ribo', 'pct_counts_ribo', 'total_counts_hb', 'pct_counts_hb', 'n_counts', 'log_counts', 'doublet_scores', 'predicted_doublets', 'doublet_info', 'cluster_0.1', 'cluster_0.2', 'cluster_0.3', 'cluster_0.4', 'cluster_0.5', 'cluster_0.6', 'cluster_0.7', 'cluster_0.8', 'cluster_0.9', 'cluster_1.0', 'cluster_1.1', 'cluster_1.2', 'cluster_1.3', 'cluster_1.4', 'cluster_1.5', 'Cell type', 'Celltype'
    var: 'gene_ids', 'feature_types', 'n_cells', 'mt', 'ribo', 'hb', 'n_cells_by_counts', 'mean_counts', 'pct_dropout_by_counts', 'total_counts', 'highly_variable', 'means', 'dispersions', 'dispersions_norm', 'mean', 'std'
    uns: 'Type_colors', 'doublet_info_colors', 'log1p', 'hvg', 'pca', 'neighbors', 'umap', 'tsne', 'leiden', 'cluster_0.2_colors', 'cluster_0.3_colors', 'cluster_0.4_colors', 'rank_genes_groups', 'cluster_0.1_colors', 'cluster_0.5_colors', 'cluster_0.6_colors', 'cluster_0.7_colors', 'cluster_0.8_colors', 'cluster_0.9_colors', 'cluster_1.0_colors', 'cluster_1.1_colors', 'cluster_1.2_colors', 'cluster_1.3_colors', 'cluster_1.4_colors', 'cluster_1.5_colors', 'Cell type_colors', 'Celltype_colors'
    obsm: 'X_pca', 'X_umap', 'X_tsne'
    varm: 'PCs'
    obsp: 'distances', 'connectivities'
In [368]:
adata.obs["Type"] = adata.obs.Type.cat.set_categories([
       'Control','Vehicle', 'Drug', 'Taxol'])
In [369]:
adata.obs["Type"] = adata.obs["Type"].astype('category')
adata.obs["Type"].cat.categories
Out[369]:
Index(['Control', 'Vehicle', 'Drug', 'Taxol'], dtype='object')
In [370]:
Sample = adata.obs['Type']
cell_type = adata.obs['Celltype']
A = adata.obs
In [372]:
%%R -i A  -o rt
rt <- A
match_celltype_levels <- c('Microglial', 'Astrocyte', 'OPC', 'ODC', 'Neuron',
       'Ependymal', 'Fibroblast', 'Endothelial', 'Dividing cell', 'Lymphocyte',
       'NK cell', 'Th17', 'B cell', 'Monocyte', 'Macrophages', 'Neutrophil',
       'Dendritic cell')
allcolour <- c('#FF0033', '#339933', '#FF6600', '#FF9999', '#990066', '#CC00CC',
       '#33FF33', '#000099', '#FFFF00', '#FF6666', '#FFCC66', '#990000',
       '#FF9933', '#0099CC', '#0066CC', '#9999FF', '#3399FF')
a <- rt
length(colnames(a))
colnames(a)[34] = c('cell.type')
a[1:10,]
mydata <- a%>%
  group_by(Type) %>%
  mutate(celltype = factor(cell.type,levels = match_celltype_levels)) %>%
  arrange(celltype)
p1 = ggplot() +
  geom_bar(data = mydata, aes(x = Type, fill = factor(celltype)), position = position_fill(reverse = TRUE)) +scale_fill_manual(values = allcolour) +
  labs(fill = "cell type",x ='Group', y = "Fraction of cells")+theme(
      panel.grid.major = element_blank(), #主网格线
      panel.grid.minor = element_blank(),
      panel.background = element_rect(fill = 'white'), #背景色
      panel.border = element_rect(fill=NA,color="black", size=0.5, linetype="solid")
        )
p1 
ggsave('./result/Ratios_barplot.pdf',p1,width = 8,height = 12,dpi = 1000)
In [373]:
cluster = {
    'Microglial cell':["Apoe",'Ctss'],
    "Astrocyte":["Aqp4","Gfap"],
    "OPC":["Olig2"],
    "ODC":["Plp1"],
    'Neuron':["Map2","Dcx"],
    'Ependymal':["Tmem212","Foxj1"],
    'Fibroblast':["Dcn","Lum"],
    'Vascular endothelial cell':["Pecam1"],
    'Macrophage':["Ms4a7","Cd74"],
    'Lymphocyte':["Cd3d"],
    'Endothelial':["Cldn5","Pecam1"],
    'Dividing cell':["Cdk1","Mki67"],
    'Lymphocyte':["Cd3d"],
    'NK cell':["Nkg7"],
    'Th17':["Rora"],
    'B cell':["Cd79a","Cd79b"],
    'Monocyte':["Ly6c2","Csf1r"],
    'Macrophages':[ "Lyz2","Cd68"],
    'Neutrophils':["S100a8","S100a9"],
    'Dendritic cell':['Ccr7',"Cd74"]  
}
marker_list = []

for k,v in cluster.items():
    if isinstance(v,list):
        marker_list.extend(v)
    else:
        marker_list.append(v)
In [374]:
sc.pl.dotplot(adata, cluster, 'Celltype', dendrogram=True,cmap='OrRd',save= 'DOTmarker1.pdf')
sc.pl.dotplot(adata,marker_list, 'Celltype',cmap='OrRd',swap_axes=True,save= 'DOTmarker2.pdf')
WARNING: dendrogram data not found (using key=dendrogram_Celltype). Running `sc.tl.dendrogram` with default parameters. For fine tuning it is recommended to run `sc.tl.dendrogram` independently.
    using 'X_pca' with n_pcs = 30
Storing dendrogram info using `.uns['dendrogram_Celltype']`
WARNING: Groups are not reordered because the `groupby` categories and the `var_group_labels` are different.
categories: Microglial, Astrocyte, OPC, etc.
var_group_labels: Microglial cell, Astrocyte, OPC, etc.
WARNING: saving figure to file result\dotplot_DOTmarker1.pdf
WARNING: saving figure to file result\dotplot_DOTmarker2.pdf
In [375]:
sc.pl.stacked_violin(adata, cluster, 'Celltype', cmap='OrRd',swap_axes=False,save='marker.pdf')
WARNING: saving figure to file result\stacked_violin_marker.pdf
In [386]:
import gc
gc.collect()
Out[386]:
145687
In [387]:
adata.write('./01.single data/scRNA_celltype.h5ad')
In [382]:
sc.settings.figdir = './result/'
sc.settings.set_figure_params(dpi=100,dpi_save=600,figsize=(5,5))
sc.pl.umap(adata,color=marker_list,ncols=4,cmap='OrRd',size=15,legend_fontsize=12,save='allmarker.pdf')
WARNING: saving figure to file result\umapallmarker.pdf
---------------------------------------------------------------------------
MemoryError                               Traceback (most recent call last)
Input In [382], in <cell line: 3>()
      1 sc.settings.figdir = './result/'
      2 sc.settings.set_figure_params(dpi=100,dpi_save=600,figsize=(5,5))
----> 3 sc.pl.umap(adata,color=marker_list,ncols=4,cmap='OrRd',size=15,legend_fontsize=12,save='allmarker.pdf')

File D:\conda\lib\site-packages\scanpy\plotting\_tools\scatterplots.py:669, in umap(adata, **kwargs)
    610 @_wraps_plot_scatter
    611 @_doc_params(
    612     adata_color_etc=doc_adata_color_etc,
   (...)
    616 )
    617 def umap(adata, **kwargs) -> Union[Axes, List[Axes], None]:
    618     """\
    619     Scatter plot in UMAP basis.
    620 
   (...)
    667     tl.umap
    668     """
--> 669     return embedding(adata, 'umap', **kwargs)

File D:\conda\lib\site-packages\scanpy\plotting\_tools\scatterplots.py:459, in embedding(adata, basis, color, gene_symbols, use_raw, sort_order, edges, edges_width, edges_color, neighbors_key, arrows, arrows_kwds, groups, components, dimensions, layer, projection, scale_factor, color_map, cmap, palette, na_color, na_in_legend, size, frameon, legend_fontsize, legend_fontweight, legend_loc, legend_fontoutline, colorbar_loc, vmax, vmin, vcenter, norm, add_outline, outline_width, outline_color, ncols, hspace, wspace, title, show, save, ax, return_fig, **kwargs)
    457     return fig
    458 axs = axs if grid else ax
--> 459 _utils.savefig_or_show(basis, show=show, save=save)
    460 if show is False:
    461     return axs

File D:\conda\lib\site-packages\scanpy\plotting\_utils.py:312, in savefig_or_show(writekey, show, dpi, ext, save)
    310 show = settings.autoshow if show is None else show
    311 if save:
--> 312     savefig(writekey, dpi=dpi, ext=ext)
    313 if show:
    314     pl.show()

File D:\conda\lib\site-packages\scanpy\plotting\_utils.py:288, in savefig(writekey, dpi, ext)
    286 # output the following msg at warning level; it's really important for the user
    287 logg.warning(f'saving figure to file {filename}')
--> 288 pl.savefig(filename, dpi=dpi, bbox_inches='tight')

File D:\conda\lib\site-packages\matplotlib\pyplot.py:958, in savefig(*args, **kwargs)
    955 @_copy_docstring_and_deprecators(Figure.savefig)
    956 def savefig(*args, **kwargs):
    957     fig = gcf()
--> 958     res = fig.savefig(*args, **kwargs)
    959     fig.canvas.draw_idle()   # need this if 'transparent=True' to reset colors
    960     return res

File D:\conda\lib\site-packages\matplotlib\figure.py:3019, in Figure.savefig(self, fname, transparent, **kwargs)
   3015     for ax in self.axes:
   3016         stack.enter_context(
   3017             ax.patch._cm_set(facecolor='none', edgecolor='none'))
-> 3019 self.canvas.print_figure(fname, **kwargs)

File D:\conda\lib\site-packages\matplotlib\backend_bases.py:2319, in FigureCanvasBase.print_figure(self, filename, dpi, facecolor, edgecolor, orientation, format, bbox_inches, pad_inches, bbox_extra_artists, backend, **kwargs)
   2315 try:
   2316     # _get_renderer may change the figure dpi (as vector formats
   2317     # force the figure dpi to 72), so we need to set it again here.
   2318     with cbook._setattr_cm(self.figure, dpi=dpi):
-> 2319         result = print_method(
   2320             filename,
   2321             facecolor=facecolor,
   2322             edgecolor=edgecolor,
   2323             orientation=orientation,
   2324             bbox_inches_restore=_bbox_inches_restore,
   2325             **kwargs)
   2326 finally:
   2327     if bbox_inches and restore_bbox:

File D:\conda\lib\site-packages\matplotlib\backend_bases.py:1648, in _check_savefig_extra_args.<locals>.wrapper(*args, **kwargs)
   1640     _api.warn_deprecated(
   1641         '3.3', name=name, removal='3.6',
   1642         message='%(name)s() got unexpected keyword argument "'
   1643                 + arg + '" which is no longer supported as of '
   1644                 '%(since)s and will become an error '
   1645                 '%(removal)s')
   1646     kwargs.pop(arg)
-> 1648 return func(*args, **kwargs)

File D:\conda\lib\site-packages\matplotlib\_api\deprecation.py:386, in delete_parameter.<locals>.wrapper(*inner_args, **inner_kwargs)
    381 @functools.wraps(func)
    382 def wrapper(*inner_args, **inner_kwargs):
    383     if len(inner_args) <= name_idx and name not in inner_kwargs:
    384         # Early return in the simple, non-deprecated case (much faster than
    385         # calling bind()).
--> 386         return func(*inner_args, **inner_kwargs)
    387     arguments = signature.bind(*inner_args, **inner_kwargs).arguments
    388     if is_varargs and arguments.get(name):

File D:\conda\lib\site-packages\matplotlib\backends\backend_pdf.py:2785, in FigureCanvasPdf.print_pdf(self, filename, dpi, bbox_inches_restore, metadata)
   2780 file.newPage(width, height)
   2781 renderer = MixedModeRenderer(
   2782     self.figure, width, height, dpi,
   2783     RendererPdf(file, dpi, height, width),
   2784     bbox_inches_restore=bbox_inches_restore)
-> 2785 self.figure.draw(renderer)
   2786 renderer.finalize()
   2787 if not isinstance(filename, PdfPages):

File D:\conda\lib\site-packages\matplotlib\artist.py:73, in _finalize_rasterization.<locals>.draw_wrapper(artist, renderer, *args, **kwargs)
     71 @wraps(draw)
     72 def draw_wrapper(artist, renderer, *args, **kwargs):
---> 73     result = draw(artist, renderer, *args, **kwargs)
     74     if renderer._rasterizing:
     75         renderer.stop_rasterizing()

File D:\conda\lib\site-packages\matplotlib\artist.py:50, in allow_rasterization.<locals>.draw_wrapper(artist, renderer)
     47     if artist.get_agg_filter() is not None:
     48         renderer.start_filter()
---> 50     return draw(artist, renderer)
     51 finally:
     52     if artist.get_agg_filter() is not None:

File D:\conda\lib\site-packages\matplotlib\figure.py:2810, in Figure.draw(self, renderer)
   2807         # ValueError can occur when resizing a window.
   2809 self.patch.draw(renderer)
-> 2810 mimage._draw_list_compositing_images(
   2811     renderer, self, artists, self.suppressComposite)
   2813 for sfig in self.subfigs:
   2814     sfig.draw(renderer)

File D:\conda\lib\site-packages\matplotlib\image.py:132, in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    130 if not_composite or not has_images:
    131     for a in artists:
--> 132         a.draw(renderer)
    133 else:
    134     # Composite any adjacent images together
    135     image_group = []

File D:\conda\lib\site-packages\matplotlib\artist.py:50, in allow_rasterization.<locals>.draw_wrapper(artist, renderer)
     47     if artist.get_agg_filter() is not None:
     48         renderer.start_filter()
---> 50     return draw(artist, renderer)
     51 finally:
     52     if artist.get_agg_filter() is not None:

File D:\conda\lib\site-packages\matplotlib\axes\_base.py:3082, in _AxesBase.draw(self, renderer)
   3079         a.draw(renderer)
   3080     renderer.stop_rasterizing()
-> 3082 mimage._draw_list_compositing_images(
   3083     renderer, self, artists, self.figure.suppressComposite)
   3085 renderer.close_group('axes')
   3086 self.stale = False

File D:\conda\lib\site-packages\matplotlib\image.py:132, in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    130 if not_composite or not has_images:
    131     for a in artists:
--> 132         a.draw(renderer)
    133 else:
    134     # Composite any adjacent images together
    135     image_group = []

File D:\conda\lib\site-packages\matplotlib\artist.py:37, in allow_rasterization.<locals>.draw_wrapper(artist, renderer)
     35 if artist.get_rasterized():
     36     if renderer._raster_depth == 0 and not renderer._rasterizing:
---> 37         renderer.start_rasterizing()
     38         renderer._rasterizing = True
     39     renderer._raster_depth += 1

File D:\conda\lib\site-packages\matplotlib\backends\backend_mixed.py:83, in MixedModeRenderer.start_rasterizing(self)
     79     r = process_figure_for_rasterizing(self.figure,
     80                                        self._bbox_inches_restore)
     81     self._bbox_inches_restore = r
---> 83 self._raster_renderer = self._raster_renderer_class(
     84     self._width*self.dpi, self._height*self.dpi, self.dpi)
     85 self._renderer = self._raster_renderer

File D:\conda\lib\site-packages\matplotlib\backends\backend_agg.py:93, in RendererAgg.__init__(self, width, height, dpi)
     91 self.width = width
     92 self.height = height
---> 93 self._renderer = _RendererAgg(int(width), int(height), dpi)
     94 self._filter_renderers = []
     96 self._update_methods()

MemoryError: In RendererAgg: Out of memory
In [ ]: